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I. INTRODUCTION 

Recent advances in numerical relativity are expand- 
ing the dynamical range covered by non-linear, three- 
dimensional simulations of general relativistic systems. 
However, we are still far from having evolutions that 
successfully simulate a binary black-hole system starting 
from the in-spiral of the holes, continuing through their 
merger and ending in the ring-down of the resulting sin- 
gle black hole. The obstacles defining this challenge are 
multiple and not all understood or even known. For ex- 
ample, it has become increasingly evident that the form 
of the Einstein equations implemented in numerical codes 
is one of the important aspects determining the behavior 
of the simulations. In principle, there are an infinite num- 
ber of ways to recast the Einstein equations as a Cauchy 
problem. The search for the formulation or family of for- 
mulations that will lead to long-lived evolutions has, in 
some situations, become a tour-de-force. In spite of the 
many "3+1 flavors" of the Einstein equations so far in- 
troduced, there has not been a formulation exhibiting a 
clear superiority. 

Besides the form of the evolution equations, there are 
other essential ingredients that have been demonstrated 
to yield improvements in the duration of the simulations. 
Gauge or coordinate conditions, methods to handle the 
singularities, non-reflective boundary conditions and dis- 
cretization methods are some examples of these elements 

m 

Another major focus in numerical relativity has been 
producing long-lasting, three-dimensional evolutions of 
the simplest case of a space-time containing black-hole 
singularities, namely that of a single, excised black hole 



in which the location of the black hole remains fixed on 
the computational grid 0- It is safe to say that this 
milestone has been achieved [3. 0- Q One should keep in 
mind, however, that some of the elements used to obtain 
simulations lasting forever cannot be extended or applied 
to binary black- hole orbits and coalescences. Specifically, 
the elements in question, as we will later discuss in detail, 
are those that hinge on keeping fixed the coordinate loca- 
tion of the black-hole singularities and thus the excision 
boundary. 

For sufficiently large separations, even in the non-linear 
regime, it should be possible to construct an approximate 
co-rotating coordinate system in which the coordinate 
drift of the black holes is avoided via an appropriate set 
of gauge or coordinate conditions ■ However, once a 
given co-rotating frame can no longer keep the location 
of the black holes fixed without introducing extreme grid 
stretching, construction of a new co-rotating frame will 
be needed. It is then likely that in this new co-rotating 
frame, the coordinate location of the black holes, and 
therefore the excision regions, will be different. Conse- 
quently, there is a natural need to develop an excision 
algorithm that is able to accommodate motion of the 
black holes. This implies handling, among other things, 
grid points that emerge from inside of the excised regions 
and become part of the computational domain. 

There are examples of simulations with moving exci- 
sion already in the literature. Using a characteristic for- 
mulation, long-lasting evolutions of a wobbling black hole 
have been obtained Sj . Other examples of simulations in- 
volving motion of the excision boundary include boosted, 
single black holes [|| and grazing collisions [ljj. -R- 0- 
cently, Yo, Baumgarte and Shapiro presented a sim- 
ple scheme for a moving excision boundary when solving 
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the scalar field equation in the background of a boosted 
black hole. 

The objective of this paper is to enhance the excision 
algorithm used in long-term stable simulations of single 
non-moving black holes to handle the case of moving 
black holes. Our approach can be viewed as an exten- 
sion of the work in Ref. [Tl| . The paper is structured 
as follows: the formulation of the Einstein equations and 
the method used in our code, called Maya, to solve these 
equations are summarized in §[n] The crux of this paper, 
namely the implementation of excision, is given in § IIIII 
Tests of the excision method for a moving black hole are 
presented in § IIVI Finally, we discuss the ramifications 
of our results and future developments in Slvl 
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II. THE MAYA CODE 

When viewed as a Cauchy problem, the Einstein equa- 
tions provide a mechanism to construct the time-history 
of the geometry of three-dimensional space-like hyper- 
surfaces. The geometry of each hypersurface is charac- 
terized by the intrinsic metric gij and extrinsic curvature 
Kij. The space-time foliation is glued together with the 
help of the lapse function a and the shift vector (3 % . In 
principle, the problem reduces to finding g^ and Kij for 
all the hypersurfaces by solving the Einstein equations 
written explicitly in terms of gij and K^, namely the 
Einstein equations in a 3+1 form. An example of these 
equations are the ADM equations 0, However, 
evolutions obtained with the ADM equations have had 
very limited success. Nakamura and Shibata 01 > an d 
later Baumgarte and Shapiro [15|. introduced a formu- 
lation (BSSN) of the Einstein equations that has clearly 
yielded remarkable improvements. Instead of g^ and 
as primary variables, the BSSN formulation introduces 
new variables g^, K, A^ and T\ The relationships 
between the BSSN and ADM variables are: 
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where Ay = K. l0 — g io K/3. Above, Eqs. Q and © im- 
ply that the conformal metric has unit determinant, 
which in turn yields the second equality in Eq. (jSJ. In 
addition to the standard constraints, Hamiltonian and 
momentum constraints, the BSSN equations require sat- 
isfying the conditions y/g = 1 and A 1 i = if one chooses 
to evolve all the tensor components of and A^ . Ex- 
plicitly enforcing ^fg = 1 is not necessary, but insuring 
that A^ remains trace-free is crucial in extending the life 
of numerical evolutions. In terms of the new variables 
(|1I5H , the Einstein evolution equations take the following 



The term in Eq. @ with superscript TF denotes the 
trace-free part of the tensor terms between brackets. 
Eq. I|10|) has been modified from the original BSSN form 
(x-terms) following Ref. 6]. This modification is crucial 
to achieve stability in single, non-moving black-hole simu- 
lations. Finally, the first terms in the r.h.s. of Eqs. H6I10|) 
involving Cp denote Lie derivatives of the scalar K and 
tensor densities gij and A^ , as well as their extension to 
quantities related to tensor densities, such as $ and T\ 
We have developed a numerical code, called Maya, to 
solve Eqs. <|6I1(J|1 . The code views this system of equa- 
tions as having the following structure: 

d t u = (3 l d t u + T(u, dp) + S{u, du, d 2 u) = p (11) 

where u = {<$>, g^, K, Aij, T 1 }. Notice that the terms 
Cpu have been explicitly separated into /3 % diU + T(u, 8(3) 
terms. In the interior of the computational domain, that 
is, excluding the outer and excision boundaries, spatial 
derivatives appearing in T and S are approximated using 
second-order, centered finite differencing. The terms of 
the form (3 l diU are commonly called "advection" terms. 
They are approximate via a second order or higher up- 
wind differencing; for details, see Ref. [l6| . 

At the outer boundary, the traditional approach for 
BSSN-based evolutions has been to assume that 
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where u represents an analytic solution. For the case 
of single black holes, this analytic solution is the exact 
solution to the Einstein equations. If an exact solution 
is not available, the flat space-time solution is used pro- 
vided the outer boundary is sufficiently far away from 
the holes. Eq. l)12ll is then implemented in a differential 
form, namely 



d t u = di(u - u ) 

r 



(u - u a ) 



(13) 



The r.h.s. of Eq. I|13|l is discretized as in the advection 
term in Eqs. l|t)ll(J|) . It was found that in order to 
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achieve stability for single, non-moving black-hole evo- 
lutions, condition l|13fl should not be used for the con- 
nection r\ For this field, one needs instead to impose 
d t T i = 0. 

Notice that in deriving Eq. I|13|) . it was assumed that 
the analytic solution is time-independent. This will not 
be the case for the moving black hole in our tests. For 
simplicity, we have chosen instead to directly use the 
Dirichlct condition u = u . This will turn out to be 
one of the two factors limiting the duration of the simu- 
lations. 

At the excision boundary, the discretization of spatial 
derivatives also requires special treatment. This consti- 
tutes the focus of the paper and is addressed in detail in 
the following section. 

Once all spatial operators in Eq. 111(1 are discretized, 
u is updated using a second-order, iterative Crank- 
Nicholson scheme [lj]- Parallelization, 10 and parameter 
manipulations in the Maya code are handled via the Cac- 
tus toolkit developed by the Albert Einstein Institute in 
Golm, Germany jig . 

Finally, several studies have shown 0, @, 0, that 
stable evolutions of single black holes are not possible 
if both the lapse function and shift vector are analytic 
functions from the exact solutions. Instead, one needs 
"driver" gauge conditions such as the 1 + log condi- 
tion for a or the T-driver for (3 l ; for a review of gauge 
conditions, see Ref. j2jj. Unfortunately, these driver- 
gauge conditions have been designed having in mind 
time-independent solutions. These gauge conditions will 
not be applicable for the excision tests of a moving black 
hole under our consideration. We are then limited to us- 
ing the exact solution to specify a and fi" 1 . Here resides 
the second factor limiting the duration of our simulations. 

III. BLACK HOLE EXCISION 

From causality considerations, one argues that it is 
completely unnecessary to evolve the system of equations 
inside the event horizon of a black hole - no physically 
meaningful information can possibly propagate out to af- 
fect the space-time outside the horizon. Doing so has 
the desirable effect of removing the black-hole singular- 
ity from the computational domain. If one chooses to ig- 
nore, or excise, the interior of a black hole, the immediate 
problem one must face is that of finding its event horizon. 
The event horizon, however, is a global object requiring 
a complete knowledge of the space-time, namely the so- 
lution to the problem itself. Unruh [2l| suggested using 
instead the apparent horizon. Apparent horizons are 
defined in terms of quantities local-in-time. That is, to 
find apparent horizons, one only needs the spatial metric 
gij and extrinsic curvature Kij of the space-like hyper- 
surfaces in the 3+1 foliation. Thus, following the history 
of apparent horizons naturally adapts to the evolution of 
the system since only information at a given instant of 
time is required. The importance of apparent horizons in 



connection with excision lies on the fact that an apparent 
horizon will always be coincident or contained within the 
event horizon. Therefore, as long as the region that is 
ignored or excised is bounded by or contained within the 
apparent horizon, this region will also be interior to the 
event horizon. 

Once the interior of an apparent horizon has been re- 
moved from the calculation, the next step to address is 
whether boundary conditions are required at the excision 
boundary. If all the fields involved in the calculation have 
at the excision boundary outgoing characteristics in the 
direction of the singularity, away from the computational 
domain, there is no need for imposing boundary condi- 
tions. In this case, the fields at the excision boundary 
are evolved in the same form as the fields in the interior 
of the computational domain. 

The problem in general relativity is that the only char- 
acteristics that one can blindly assume to be outgoing at 
the excision boundary are those of physical modes since 
only those modes are causally constrained. For other 
modes, such as gauge modes, the characteristics depend 
on the particular structure of the evolution equations 
as well as the set of gauge conditions imposed. Mani- 
festly hyperbolic formulations of the Einstein equations 
facilitate determining the characteristics of all the fields, 
thus allowing the identification of those fields requiring 
conditions at the excision boundary. Here is one of the 
main reasons for the popularity of explicitly hyperbolic 
formulations of the Einstein equations p3| . This would 
suggest that numerical simulations involving excision are 
only possible using a hyperbolic formulation. However, 
studies by several groups using the BSSN formulation, a 
formulation that is not explicitly hyperbolic, have shown 
[H miHt that this is not the case. BSSN-based codes, 
assuming all fields having outgoing characteristics at the 
excision boundary, have produced results similar to those 
based on hyperbolic formulations. 

In summary, we work under the assumption that there 
is no need for boundary conditions at the excision bound- 
ary. The only task is to design a discretization of the 
evolution equations that is appropriate at the excision 
boundary where centered finite differencing is no longer 
feasible. 

A possible approach to handle spatial derivatives at 
the excision boundary is to modify all finite-difference 
stencils at this boundary so as to avoid using grid-points 
from within the excised region (e.g. one-sided differ- 
ences), while still keeping the same order of accuracy. 
This typically involves getting lost in a veritable jun- 
gle of logic, especially when the excision region is that 
of the "Lego sphere" , a representation of a sphere us- 
ing only the points on a fixed Cartesian grid. Each type 
of corner must be treated with its own finite difference 
stencil. This approach was followed in early simulations 
0, ^3 ■ To circumvent the complications of this stencil- 
modification method, Alcubierre and Brugmann £| in- 
troduced a very simple excision method which is ideally 
suited for stationary space-times such as those involved 
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with single black- hole evolutions. In simulations of 
single black holes that did not assume any symmetries 
had limited duration. It has been recently demonstrated 
0, however, that the factor restricting the duration of 
these simulations was not the excision method used, but 
a term in the evolution equation of the connection T l . 

Anticipating that for more general cases, e.g. binary 
black holes, the excision boundary and the r.h.s. in 
Eqs. H6I10|I w hl have strong temporal dependence, we 
present a method for handling the excision boundary that 
generalizes to non-stationary space-times. Our approach 
extends the method used by Yo, Baumgarte and Shapiro 
in their study of the dynamics of a scalar field in the 
background of a single, boosted black hole . We have 
carried out tests to show that our method, when applied 
to the case of a non-moving black hole, reproduces the 
observed long-term stability observed in Ref . Q . 

Once an apparent horizon is found, the first step is to 
choose the shape and size of the excision region. A brief 
comment is warranted regarding the choice of the shape 
of the excision region. Historically, a sphere was used. 
The rationale for this is that the apparent horizon will 
tend to be spherical or ellipsoidal in all but extremely dy- 
namical space-times. As mentioned before, on the Carte- 
sian grids typically used in three-dimensional numerics, 
one achieves only a crude approximation of the sphere, 
colloquially known as a "Lego sphere". Alcubierre and 
Brugmann simplified matters by using a cube as an 
excision shape. We have experimented with a variety of 
excision shapes, all possessing at the very least octahe- 
dral symmetry (because of the frequent use of octant, 
quadrant, and equatorial symmetry boundary conditions 
in performing three-dimensional simulations) , and in the 
end, have settled on the sphere, the cube, and the cuboc- 
tahedron. These shapes are shown in Fig. ^ 




FIG. 1: The "idealized" excision shapes one may wish to 
use for excision - a cube, a cuboctahedron, and a sphere, 
along with their representations on a Cartesian grid. The 
commonly-seen "Lego sphere" at the bottom right is the 
largest of the three shapes in terms of the volume it encloses. 

To facilitate the excision method, we adopt the com- 
mon practice of carrying an extra grid function, called the 
"mask" , to indicate the state of a grid-point in the com- 
putational domain. The typical states of a grid-point are: 
excised, excision boundary, interior, and outer boundary. 
Fig. [21 depicts a schematic representation of the mask 



variable. The large, dark circle is the location of the ap- 
parent horizon. The excision region is a "Lego sphere" 
represented by the thick line. The filled circles are the 
excision boundary points, while the empty circles are the 
excised points. We also use the mask function to label 
the outer boundary points, but this is irrelevant in the 
present discussion on excision. 




FIG. 2: Schematic of the mask function for the case of spher- 
ical excision approximated on a Cartesian grid. The large, 
dark circle represents the location of the apparent horizon, 
while the light circle is the spherical excision region which is 
approximated by the dark line. The small, filled circles are 
the points on the excision boundary. 

Given u n at the interior, excision boundary and outer 
boundary points, the r.h.s. p of the evolution equations 
(|ll|l for the interior grid-points are constructed and used 



to update u n 



at these points. Before any subse- 



quent updating takes place, including intermediate iter- 
ated Crank-Nicholson steps, updating of u n at the exci- 
sion boundary is required. 

Our approach for updating data at the excision bound- 
ary is designed to avoid altering the finite difference sten- 
cils used for the interior points. This is crucial for pre- 
serving a simple code. It consists of obtaining the data 
at the excision boundary via polynomial extrapolation. 
In order to simplify matters even further, we perform 
one-dimensional extrapolations along normal directions 
to the excision boundary. Because in general the ex- 
cision boundary is only a discretized approximation to 
the geometric shape desired, the normal to the excision 
boundary is not aligned with the grid-points in the com- 
putational domain. We select the direction of extrapo- 
lation to be as close as possible to the direction of the 
normal for which all data involved on the extrapolation 
lie on grid-points. 

We consider two approaches to construct data at the 
excision boundary. One approach involves the extrapola- 
tion of p from the interior points to the excision boundary 
(rhs-extrapolation method) . The other approach is based 
on extrapolation of the solution u n+1 (sol-extrapolation 
method). Fig. [3] shows an illustration of both methods. 
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Excised points are labeled with an X, the excision bound- 
ary points with a box and interior points with a solid dot. 
As mentioned before, for each interior point we can use 
centered finite differencing to calculate p. That is, P2 
through pq in Fig. [3| are known, but p\ is unknown. The 
value of pi is, in principle, required to obtain The 
rhs-extrapolation method finds the value pi using the 
values of p2 through ps, depending on the order of the 
extrapolation. 

Given that for the interior points, the truncation error 
in the discretizations to obtain p is second-order, the rhs- 
extrapolation must be such that the overall truncation 
error is preserved. This implies rhs-extrapolations of at 
least second-order. That is, 



pi=2p 2 - p 3 + 0(h 2 ). 



(14) 



with h the grid spacing, ft is important to notice that 
the simple excision method first used by and later by 
consists of first-order extrapolation, namely p\ = p2- 
This choice yields truncation errors of 0(h). Our numer- 
ical experiments with a moving black hole using zero- 
order extrapolations quickly became unstable, in agree- 
ment with similar results in A possible reason 
why zero-order extrapolation works for non-moving sin- 
gle black holes is that, from Eq. I|llfl . the p\ — P2 condi- 
tion is equivalent in the continuum to 



= d r p = d r d t u = d t d r u , 



(15) 



which is consistent with the time-independent nature of 
the solution. 

The sol-extrapolation method uses the values u^ +1 
through to obtain as depicted by the solid 

black arrow. In order to preserve second-order accuracy 
of the solution in this case, third-order extrapolations are 
needed. This is because the extrapolated value it™ +1 will 



depicts this process. The movement is indicated by the 
shift in the apparent horizon (dark circles), which causes 
the excision region to shift in order to remain centered on 
the apparent horizon. The dotted circles are the original 
apparent horizon and excision region. The filled circles 
are the new excision boundary points that were previ- 
ously excised points. While nothing needs to be done 
with interior points, we need to populate the previously 
excised points with values. This is done following the 
sol-extrapolation method previously described. 




be used in computing P2 for the updating u% 



Tt + l 



U 



n+2 



FIG. 4: This figure is similar to Fig. [5] however, now the 
excision boundary has moved, shown here by the shift in the 
apparent horizon. The excision region must shift to remain 
centered within the apparent horizon, thus causing a need to 
populate previously excised points, represented here as filled 
circles. 



IV. DANCING BLACK HOLE 



Mill 




FIG. 3: In this schematic illustration of extrapolation onto 
the excision boundary, an X represents an excised point, a 
box represents an excision boundary point, and a solid dot 
represents an interior point. 

As the coordinate location of a black hole changes, so 
does the position of the excision region. We adjust the 
numerical evolution such that the excision region does 
not move by more than a grid-point. The consequences 
of this change of location are that points previously la- 
beled excision boundary are now interior and points pre- 
viously labeled excised are now excision boundary. Fig.0] 



In this section, we present numerical experiments that 
demonstrate the ability of the Maya code to handle mov- 
ing excision boundaries. Because the algorithms and data 
structures presented here are general, they can provide 
a preliminary treatment for non-trivial time-dependent 
systems such as black- hole binaries. 

In order to test our excision algorithm, we introduce 
the simplest time-dependency into the system that does 
not change the form of the metric. We consider a single, 
non-rotating black hole in ingoing Eddington-Finkelstein 
coordinates and apply a coordinate transformation of the 
form: 



where 



d l = 
d l = 



d 1 . 



A sin [ ,0,0 or 
S.,,s ( ,,4 sin [^Pj ,0 



(16) 

(17) 
(18) 



6 



Notice that this is a coordinate transformation of the spa- 
tial coordinates only. The above transformation induces 
the following change on the shift vector: 



where 



b l = 
V = 



/3 l -> (3 l + V 



[ ^ ) ■<».() 



or 
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(19) 

(20) 
(21) 



The functional form of a, gij and Kij remains the same. 
The only effect of this coordinate transformation is to 
change the coordinate location of the black hole. For 
the case (|17f) . the black hole bounces back and forth 
with an amplitude A and maximum coordinate veloc- 
ity V (bouncing- hole). Similarly, for the case l|18l) . the 
black hole moves in a circle of radius A and coordinate 
velocity V (circling- hole) . 

The size of the excision region depends on the shape 
chosen. In the case of excising a spherical region, the 
excision radius was set to 1.5M with M the mass of the 
black hole and 2M its radius, while the cubical excision 
region is inscribed to fit inside a sphere of radius 1.5M. 
That is, we are excising a larger region of the computa- 
tional domain when we use a sphere. We have empirically 
found that sol-extrapolation performs a little better than 
rhs-extrapolation for the tests under consideration. Wc 
use third-order (0(h 3 )) sol-extrapolation for the excision 
boundary and fourth-order for repopulating. The sim- 
ulations were carried out with grid-spacings h = 0.2M 
and imposing equatorial symmetry about the z-axis. The 
computational domain has a total length of 20M in the 
x and y directions and 7M in the z direction. We used 
values of A = 1.0M and V — 0.5 for the simulations 
presented here. 

The driver gauge conditions commonly used to pro- 
duce long-term stability in single, non-moving black-hole 
evolutions are the 1 + log slicing and T-driver; that is, 



= \d t r, 



(22) 
(23) 



respectively. Unfortunately, these gauge condition do not 
work for our moving black-hole tests. The black-hole so- 
lution that one obtains under the coordinate transforma- 
tion H16| l are not compatible with these conditions. New 
or modified driver gauge conditions would have to be con- 
structed. These new driver conditions would likely be 
dependent on the particular nature of our moving black- 
hole solution, thus not generalizable to astrophysically 
interesting cases such as black-hole binaries. Because the 
focus of the work in this paper is on excision, we have 
chosen to set a and (3 l to the analytic values given by 
the exact solution. It is well known that using analytic 
lapse and shift gauge conditions yields evolutions of sin- 
gle, non- moving black hole lasting < lOOAf . As we will 



see next, the limits in duration of our moving black-hole 
simulations are consistent with the stationary cases. 

As mentioned before, a similar problem arises with the 
boundary condition (|13fl . This condition assumes that 
the analytic solution is time- independent. For the case 
of time-dependent solutions, the condition becomes 

d t u = d tUo - -di(u - u ) + {U - Uo) . (24) 
r r 

The modified condition (|24J) should in principle work for 
our time-dependent case. However, one needs to remem- 
ber that the outer boundary condition for the connection 
T l that yields long-term stability for non-moving holes 
was not ffTT^f) . but d t T l = 0. We were unable to obtain 
a condition applicable to our time-dependent case that 
would lead to significant improvements in the duration of 
the evolutions. Therefore, for simplicity, we have chosen 
instead to directly use for all fields the analytic solution 
as the boundary condition. 

Fig. is a xt-p\ot of the variable $ for the bouncing- 
hole simulation. The sinusoidal flat region is a x-cross 
section of the world-tube of the excision boundary. Fig. El 
is also a xt-plot but in this case for the normalized Hamil- 
tonian constraint. The reflections at the outer boundary 
are clear from this figure. These effects become stronger 
with time and eventually are one of the main causes of 
the simulation terminating. The other contributing fac- 
tors are the analytic gauge conditions. The circling case 
is depicted in Fig. The figure consists of stacked time 
slices of the variable $. The time slices are separated 
by 15M increments in time. For clarity, the domain dis- 
played only represents ±5M in the x and y dimensions. 
The "orbit" of the black hole is Air, so each slice repre- 
sents an evolution of approximately one and a quarter 
orbits from the slice preceding it. Animations of both 
bouncing and circling holes can be found at [2(J. 

For comparison of evolutions, wc monitor the Hamil- 
tonian constraint, the area of the apparent horizon, and 
a "mass" function M$ obtained from 



M$ = -e 
2 



12<I> 



1 . 



(25) 



For the exact solution, M$ is equal to the mass of the 
black hole. In all the following figures (Tigs. I9TTT)) . the 
top panels depict results using a spherical excision region 
and the bottom cubical. We also show results from a 
stationary black-hole case with a solid line. This is our 
reference case as we do not expect the moving cases to 
supersede the stationary case. The dotted lines represent 
the bouncing-hole and the dashed the circling-hole. 

Fig. |5| contains plots of the L2-norm of the Hamilto- 
nian constraint versus time for each case. It is evident 
from this figure that the simulation is unstable even in 
the case of a stationary black hole. The instability is 
likely due to the gauge and boundary conditions. When 
we turned on the driver conditions IL'l'l) and l|23|) , as well 
as the boundary condition (|13fl . long-term stable runs 
for a stationary hole were obtained, as one can see from 
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FIG. 5: Space-time plot of the BSSN variable $ for a 
bouncing-hole. The plot shows the evolution in the interval 
< t < 25M. Time runs along the sinusoidal canyon. 




FIG. 6: Space-time plot of the normalized Hamiltonian con- 
straint for bouncing-hole. The plot shows the evolution in the 
interval < t < 25M. Time runs along the sinusoidal canyon. 



1 




FIG. 7: Stacked time slices of the BSSN variable <E> for a 
circling-hole, with time proceeding vertically upward, and the 
time slices separated by 15M. The domain displayed only 
represents ±5M in the x and y dimensions, for greater clarity. 
The "orbit" of the black hole is 4tt, so each slice represents 
an evolution of approximately one and a quarter orbits from 
the slice preceding it. 

the jagged behavior is because the L2-norm of M$ was 
done with data that did not include excision boundary 
data. The data were taken from within a shell centered 
at the location of the black hole a distance 0.5 M from 
the excision boundary. The jaggedness in the Hamilto- 
nian constraint plots is due to the artificial jump of values 
because the center of the excision region is forced to be 
at a grid point, thus producing discontinuous motion of 
the excision region. 



Fig. [SJ similar to those in Fig. 1101 shows the area of 
the apparent horizon versus time; the apparent horizon 
tracker is described in |2?| . The apparent horizon tracker 
fails to give reasonable results at approximately t = 90M 
because the coordinate distortions are such that the hori- 
zon intersects the excision region at one or more points. 
Finally, Fig.llllshows plots of the L2-norm of M$ versus 
time. The reason why the L2-norm of Af$ does not show 



V. CONCLUSIONS 

Ref |4j described a simple excision algorithm that ex- 
tended the lifetimes of time-independent, single black- 
hole evolutions. This work was expanded in Ref. [Tlj| . 
showing that higher order extrapolation is necessary for 
stability in studying the dynamics of a scalar field in the 
background of a boosted Kerr-Schild background. We 
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FIG. 8: L2-norm of the time derivative of K for a single, sta- 
tionary black hole with driver conditions 1221231 and bound- 
ary condition 1131 . 



id 
a; 
u 

< 



51 



50 



49 



60 - 



55 



50 



1 cubical excision 




FIG. 10: Same as in Fig. but for the area of the black 
hole's apparent horizon. For reference, the area of a single 
black hole is 16 n M 2 = 50.2655 M. 
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FIG. 9: L2-norm of the Hamiltonian constraint. The top 
panel shows results when a cubical excision is used and the 
bottom panel corresponds to spherical excision. In each of the 
panels the solid line depicts data from a stationary black hole, 
the dotted line from a bouncing black-hole and the dashed line 
from a circling black-hole. 
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FIG. 11: Same as in Fig. El but for the L2-norm of the mass 
function M*. 



have generalized and applied the excision algorithm pre- 
sented in ^l| to the case of a single black hole moving 
its coordinate location throughout the computational do- 
main. We showed that this algorithm is generalizable to 
non-stationary black-hole space-times. We also find that, 
despite the fact that a cubical excision region is all that 
is necessary for long-lived stationary black-hole runs, a 
spherical excision boundary is viable. A spherical shape 
may be more desirable in dynamic simulations because 
it allows the excision region to fit more easily within the 



apparent horizon. It is important to emphasize that the 
algorithm in this work for excising the black hole as it 
moves through the computational domain behaves in a 
similar way to excising a stationary black hole. This 
suggest that as gauge and boundary conditions suitable 
for time-dependent solutions are developed, the excision 
method discussed in this work can be directly imple- 
mented without requiring further or with minimal mod- 
ifications. 
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